function [drift, bellows_volt_drifted] = drift_correction(patient);
    %|------------------------------------------------------------------------|
    %|Calculation of linear drift in the bellows voltage, based on the        |
    %|minimization of the model error, using fminsearch to vary drift.        |
    %|  Drift variable should be used to correct bellows voltage as follows;  |
    %|
    %|      bellows_volt_drifted = bellows_volt + drift*scan_times_zerod;     |
    %|
    %|  where   bellows_volt        =       originial bellows signal          |
    %|          scan_times_zerod    =       time of bellows_volt (in secs),   |
    %|                                      starting at 0s.
    %|  Runs the model on a sub-sampling of the image data, defined by;
    %|      step=10;
    %|      ii = round(dim(1)/2);  % usually just run on a single coronal slice,
    %|                             % for speed, change to
    %|                             % ii = 1:step:dim(1) for completness
    %|      jj = 1:step:dim(2);
    %|      kk = 1:step:dim(3);
    %|
    %| drift.txt (variable) is saved to model_folder.
    %------------------------------------------------------------------------
    %|      Dependancies;                                                     |
    %|            choose_patient (finds data folders)                         |
    %------------------------------------------------------------------------
    %   This file is part of the
    %   5D-Novel4DCT Toolbox  ("Novel4DCT-Toolbox")
    %   DH Thomas, Ph.D
    %   University of California, Los Angeles
    %   Contact: mailto:dhthomas@mednet.ucla.edu
    %------------------------------------------------------------------------
    % $Author: DHThomas $	$Date: 2014/04/01 10:23:59 $	$Revision: 0.1 $
    
    %Check if drift.txt exists;
    if exist([patient.model_folder '/drift.txt'])>0;
        
        drift = load([patient.model_folder '/drift.txt']);
        bellows_volt_drifted=patient.bellows_volt+drift*patient.scan_times_zerod;
        %if no, run drift correction code;
    else
        
        step=10;    %chose sub-sampling of image data to run the model, for speed.
        ii=round(patient.dim(1)/2); %1:step:patient.dim(1);%    % usually just run on a single coronal slice, for speed;  
        jj=1:step:patient.dim(2);       % solve the model every *step* voxels.
        kk=1:step:patient.dim(3);       % solve the model every *step* voxels.
        
        % define optimisation parameters;
        optionsforme = optimset('Algorithm','interior-point','MaxFunEvals',10000,'display','iter','UseParallel','never');%,'FunValCheck','on');%,'UseParallel','always');%,'MaxFunEvals',2500);
        initparam =-ones(1,1)*1E-3; %initial parameters;
        
        % using fminsearch to optimise the drift value, runs the model codel
        % from the following script;
        tic
        %     [patient.drift, temp2]=fminsearch(@optimise_drift,initparam,optionsforme,patient.folder,patient.model_folder,patient.folder_static,patient.static_mask,patient.scans,patient.run_scans,patient.bellows_volt,patient.flow,scan_times,dim,ii,jj,kk,step,leave_out_num,0,folder_output, ref);
        [drift, temp2]=fminsearch(@optimise_drift_toolbox,initparam,optionsforme,patient,ii,jj,kk);
        
        toc
        
        bellows_volt_drifted=patient.bellows_volt+drift*patient.scan_times_zerod;
        save([patient.model_folder '/drift.txt'], 'drift' ,'-ascii')
        
    end
    display(sprintf('Drift Correction = %.2f [mV/s]',drift*1E3))
    
    
    
    
    
    
